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026803 (2011)]. It is confirmed that the number of the steps corresponds to the number of stream 
lines of electrons through the constriction. We verify that density fiuctuation plays a major role in 
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I. INTRODUCTION 

Electric conduction in confined geometries has been an 
important subject of study to characterize properties in 
the mesoscopic scale [Tl], and its study has a wide pos- 
sibility of applications to electronics 0, Q . For rarefied 
gases under one constriction, such as a point contact or 
a circular orifice, the transport shows the quantized con- 
ductance [l|, II, 111 , Sharvin resistance Maxwell 
resistance [9|, Il0| . depending on the constriction width 
and the electron mean free path ^M^- 

On the other 

hand, for dense liquids in which the exclusion effects be- 
tween particles are important, the formation of layers 
[is, 19], pinning and depinning ^2Q] . and anisotropic and 
nonuniform mobility [2ll423j are observed in microchan- 
nels. However, the geometrical effects to conduction in 
such strongly correlated particle systems have not been 
fully understood yet. 

Recently, Rees et al. have found the existence of 
steplike conductance in the conduction of classical two- 
dimensional (2D) electrons on liquid ''He (Refs. [l^l and 
[25j ) using a device with point-contact geometry [26| . 
They suggest that the origins of the conductance are 
not quantum effects, as in the conductance quantization 
but are the effects of strong correlations in classical 
electron systems. We find it necessary to justify their 
suggestion and to clarify the mechanisms of the dynam- 
ics behind the steplike conduction of electrons on liquid 
^He. The main goal of this paper, thus, is to understand 
the mechanisms of the steplike conduction in classical 2D 
electron systems. 

In this paper, we investigate the electric conduction of 
the classical 2D electron system through a narrow con- 
striction in a microchannel using a molecular dynamics 
(MD) simulation with two particle baths under a Nose- 
Hoover thermostat. The electrons in the MD simulations 
are confined in a point-contact-like shape and interact 
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with each other, in terms of electric potential derived 
from the Poisson equation under a boundary condition 
given so as to imitate the device in Ref. [1^ (see Fig. [iJ. 
Our model and method stand on the presumption that 
the many-body effects in the confined geometry are es- 
sential for the steplike electric conduction (see Sec. Ill Al 
and Appendix \^ . We calculate conductance as a func- 
tion of the confinement strength of the constriction. To 
confirm the suggestion by Rees et al. ,26;] and to develop 
the understanding of the steplike conduction, we investi- 
gate the static and dynamical properties of electrons near 
the constriction from the spatial distribution of electron 
density, electrostatic potential, and potential fluctuation. 

The organization of the paper is as follows. In Sec. [Hi 
our model, methods, and simulation details are provided. 
The conductance and the other quantities calculated in 
our simulations are presented in Sec. IIIII in which the 
suggestion by Rees et al. [2^ is also verified. Based on 
the results in Sec. 1111) we discuss the mechanisms causing 
the observed conductance in Sec. IIV Al In Sec. IIVBI we 
compare our results with the observed conductance in the 
experiment. We also present the result of the Langevin 
dynamics (LD) simulation method which is the molecu- 
lar dynamics using the Langevin thermostat in Sec. IIV CI 
The results and the discussions are summarized in Sec. 
IVl In appendices, we present some detailed descriptions 
of our model and method. In Appendix \^ we estimate 
the resistance from the electron-helium vapor atom scat- 
tering and from the electron-ripplon scattering. In Ap- 
pendix |Bl we give the details of the method to calculate 
physical quantities. In Appendix [Cl we briefly explain 
the method of our LD simulation. 



II. MODEL, AND NUMERICAL METHODS 

In this section, we introduce our model with a demon- 
stration to justify our treatment for the electrons in the 
device in Ref. ^2^] and explain the method of the MD, 
including how to provide the simulation setup and to cal- 
culate observed quantities. 
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A. Electron states and transport properties in the 
electron system over liquid ''He on metal electrodes 

In this section, we briefly explain the electron state and 
the transport process of electrons in the device in Ref. 
[26j . The system of electrons over liquid "^He on metal 
electrodes studied by Rees et al. [26*1 is an ideal system 
to investigate strongly correlated classical 2D electron 
systems |23l - l28| . The electrons are confined on a plane 
at the height zq from the liquid "^He due to a poten- 
tial barrier at the liquid surface, an interaction between 
the electrons and polarized liquid ''He, and charges in- 
duced in the metal electrodes [1^. We also apply a 
holding electric field E± normal to the liquid surface. 
For surface density Ug ^ 10^ — lO^cm^^, temperature 
T ~ 1 K, and liquid thickness zhc ^ 1 A*ni, the verti- 
cal motion of each electron can be decoupled with that 
in the parallel direction to the surface. For the verti- 
cal direction, the electrons occupy the ground state of 
the quantized electron states, whereas the parallel mo- 
tion is classical. These treatments of electron motion on 
the plane can be justified because the interelectron sep- 
aration, Vs — 2/^TTns ^ 10^' /im, is much larger than 
the thermal de Broglie wavelength. Ad ~ 10~^ /im, and 
the energy gap between the ground and the first excited 
states, A ~ 19 K, is sufficiently larger than the temper- 
ature in the experiment. Here, A is estimated by the 
variational method |30] under E± ~ 6.5 x 10^ V/cm in a 
very similar device [28] to that used in Ref. [2^. 

For the bulk 2D electrons, it is known that the mo- 
bility on the surface is little affected by the scatter- 
ings by roughness on the interface liquid substrate for 
zhc ~ 1 Min [31,-,34|. but is dominated by the scatterings 
between electrons and helium gas atoms for T > 1 K 
or between electrons and ripplons at lower temperature 
[23, 24, 29, 35J. The correlation effects between electrons 
on the transport have been discussed by the kinetic equa- 
tion method under the complete control approximation 
[36l - l40j , the force- balance method [IJ, Hi], \^ , and more 
sophisticated theory '13, 'is', "i?] . 

In Appendix |^ we briefly estimate the resistance for 
the transport in the channel geometry as in Ref. p6| . 
The result supports that the many-body effects in the 
confined geometry are dominant for the steplike electric 
conduction. 



B. Interaction and confining potentials 

In this section, we present the forms of interaction po- 
tential between electrons and confining potential for our 
MD. The electric potential energy for the i— th electron 
on liquid ^He is given by [isj 

</)(r„0 = ^ </>! [njit)] + 0c [r^it)] , (1) 

where ri(<) = [xi{t),yi{t)] is the position of the i— th 
electron a.t z = Zh at time t, rij{t) = \ri{t) — rj{t)\, (j)i{r) 



is the interaction potential energy between two surface 
electrons, and (/)c(r) is the confining potential energy. We 
solve the 3D Poisson equation for (p in the semi-infinite 
domain z > 0, and presume the electrons to be confined 
at the height z = Zh from the plane z = 0. Here, we set 
Zfi as the thickness of liquid "^He, Zh = 1.5 /xm, which 
follows the experimental setup [2a] (Fig. [T]), where the 
average distance zq ~ 10~^ /im of the hovering electrons 
from the liquid surface is disregarded [13] . Therefore, the 
potential is obtained from an analytic solution at z = Zh 
of the Poisson equation. Then, (r) is given by [H, \^ 
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where e is the elementary electric charge. The right-hand 
side in Eq. ([5]) consists of the bare Coulomb interac- 
tion and the dominant screening effects between surface 
electrons which represents the contribution of the image 
charge induced in the metal electrodes. Equation ^ is 
obtained under the same order approximation to the po- 
tential energy in the previous studies [1^ [s^j . 

The boundary conditions we impose on the Poisson 
equation are 
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where Sq, Sg, and the outside of U S'g at z = 0, 
respectively, represent the reservoir, the split gate, and 
the guard electrodes [l^ [see Fig. [HJa)]. Then, (/)c(r) is 
represented as 
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FIG. 1. (Color online) A schematic view of our simulation 
setup of a classical 2D electron system confined on the plane 
a,t z — Zh- On the boundary plane z = 0, the yellow and red 
regions, respectively, correspond to the reservoir and the split- 
gate electrodes in the experiment [2^. The system consists 
of a left and a right reservoirs, and a channel of length 2D, 
where the reservoirs are connected with the channel. 
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where Vq and Vg are the voltages of a reservoir and a 
split-gate, respectively. It should be noted that each in- 
tegration in Eq. ^ can be performed exactly if we as- 
sume that the integration range consists of rectangles and 
triangles, as in Fig. [2ja). 

We also adopt Vq = 0.38 V and Vg of the range from 
-0.05 V to 0.38 V, which imitates the device in Ref. 
When Vg is set to Vb, 0c works as the confinement with- 
out the point contact, and when Vg is set to a voltage 
lower than Vq, <pc additionally generates the point con- 
tact due to a voltage induced between 5*0 and ^g, as seen 
in Figs. Mh) and (c). 



(a) 




(b) 
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FIG. 2. (Color online) (a) The boundary condition at 2; = 0: 
So and So, respectively, are the yellow and red regions. The 
electric potential in So and the outside of the colored regions 
are Vb and 0, respectively. The gate voltage Vg is imposed to 
Sq. (b) The confining potential energy 0c for Vb — 0.38 V, 
Vg = 0.38 V, and Zh = 1.5 [s,ra in Eq. ©. (c) The contour 
plot of 00 in the channel, which is the region oi \y\ < D = 
15 /im. 



lows: This method is the established one to reproduce the 
precise equilibrium state, and the equation of motion for 
electrons can keep the local time-reversal symmetry. In 
CMD, we introduce a fractional particle characterized by 
an extended number variable (ENV) for a particle bath, 
where the integer part of the ENV denotes the number 
of particles, and the fractional part represents the value 
of the fractional particle. The ENV couples with the sys- 
tem through the fractional particle interacting with the 
rest of the system. 

In our simulation, the system consists of a left and 
a right reservoir, and a channel of length ID with D = 
15 /xm, where the reservoirs are connected with the chan- 
nel (see Figs. [T]and[2]). The left (right) reservoir has an 
electrochemical potential (mr) which can be divided 
into two parts as. 



(4) 



with 7 = L, R, where /x!^ is the intrinsic part of the ideal 
chemical potential [ssj . and fj,l^, which is the control pa- 
rameter in our MD, is the sum of the excess chemical 
potential and the confining potential energy. The left 
(right) reservoir consists of one fractional particle and 
temporally variational Ni^ (^r) electrons in Sl (5'r), 
where Sl (-S'r) is the region y < —H [y > H) with 
= 15 /zm. 

The equation of motion for the electrons in our MD is 
given by 
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where 



and 



'L(R) 



are, respectively, the 



componentia — x,y]of the position of the z— th 
electron and the fractional particle belonging to the left 
(right) reservoir, rjL(R) = jr^ - rL(R)|, i^l (j^r) is the 
ENV of the left (right) reservoir, m is the electron mass, 
and 0{x) is a step function, i.e., 9{x) ~ 1 for x > 1, and 
6{x) = otherwise. Here, the "friction" coefficient C is 
adjusted according to the following equation [sot : 
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C. Constant temperature and chemical potential 
molecular dynamics 

Our MD is built on a hybrid scheme of a constant tem- 
perature MD with a Nose-Hoover thermostat[48-50] and 
a constant chemical potential MD (CMD) [5ll454j . The 
reason we adopt the Nose-Hoover thermostat is as fol- 



where is the "mass parameter" of N is the tem- 
porally variational total number of electrons, kg is the 
Boltzmann constant, and Tk is the expected kinetic tem- 
perature. The time evolutions of the fractional particle 
coordinate r" and the ENV are respectively given by 
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where Qy is the mass parameter of v^. 

The temporal variations of and Nyi are governed 
by the following protocol [HJ]- When — iV^ becomes 
zero, we delete the fractional particle, and replace one 
electron in the reservoir by a new fractional particle sat- 
isfying v'^'^ — i)!^™, when the coordinate and the velocity 
for the deleted fractional particle are discarded, and those 
of the replaced electron are given over to the new frac- 
tional particle. When — becomes one, we convert 
the fractional particle into an electron, and insert a new 
fractional particle at the position where the fractional 
particle satisfies 1)°^^ — j>"°'" and the potential energy 



also satisfies a local minimum condition i^^ 



0, when 



the coordinate and the velocity of the converted frac- 
tional particle are given over to the new electron, and 
the velocity of the new fractional particle is set to zero. 



The condition i)'. 



.old 



ensures the temporal conti- 



nuity of iy-y, and i^^ ~ works so as not to change the 
average velocity of the system due to the insertions of a 
fractional particle. Moreover, the inserting place of the 
fractional particle in the left (right) reservoir is selected 
from the squares created by dividing area Al (Ar) in 
5l {Sh) into 0.01 /xm square mesh, based on the condi- 
tion described above. The behavior of electrons in the 
channel is not disturbed by the fractional particles, be- 
cause if the fractional particles try to enter the channel, 
i)j increases, and thus the fractional particle is converted 
into an electron according to the above protocol for the 
particle conversion. 



D. Simulation setup 

Throughout our MD simulations the initial total num- 
ber of electrons Ni is 1284, and the initial positions 
of electrons are located in the ground-state configura- 
tion which forms a classical 2D Wigner crystal deformed 
by the confinement at each Vq. The initial velocities 
of electrons are randomly assigned from the Maxwell- 
Boltzmann distribution at Tk = 1.2 K corresponding to 
the experimental setup [2^. The initial positions of the 
fractional particles are determined based on the protocol 
in the case of i^-y = 0, and their initial velocities are set 
to zero, = 0. The extension variables are initially set 
to: — Nj + 0.5, i>j — 0, — 0, and C = 0. We choose 
as: Qc = ksTKNir^ = 6.02 x lO^^o meV-s^, where 
tn = 2.13 X 10^^^ is the Nose-Hoover relaxation time. 
The value of tn is to make the thermostat work effectively 
[56| . and is close to the characteristic period of short- 
wavelength vibrations of the bulk 2D Wigner crystal [s^l , 
Tg = 8.52 X 10^^^ s. It is believed that the temporal varia- 
tion of temperature does not affect long- wavelength con- 
ductivity directly because the temperature fluctuation in- 
duces only short-wavelength fluctuations in the momen- 
tum of electrons. We also carry out the MD simulations 



with another mass parameter Q'^ = IGQ^, and then ob- 
tain quantitatively similar results for all quantities cal- 
culated in the following. Although it is known that Q,y 
affects the cycle of fluctuation in the number of particles, 
we only use one parameter Q^, — Eotg, where Eq — 1.44 
meV and to ~ 6.28 x 10^^^ s are our MD units of energy 
and time, respectively. Because of the restriction of our 
computer resources, we have to limit the area Al (Ar ) to 

being the region -ymax < y < ~ym\n (Vmin <V < ^max) 

with ?/niin = 31.7 /im and J/max = 32.2 jjum in \x\ < 7 fim, 
where — ?/min (?/min) corrcsponds to the average position 
of the electron present at the furthest left (right) in all of 
the electrons in the left (right) reservoir. The above ini- 
tial conditions and the parameters are used throughout 
our simulations. 

We fix and fi}^ as fi}^ = 49.16 meV and fi^ — 48.96 
meV, where we set the minimum of the confining po- 
tential energy to zero. Thus the chemical-potential dif- 
ference Afi^ = /J.L ^ /^R induces a direct current (dc) 
between the reservoirs. Although Rees et al. [2^ mea- 
sured an alternative current (ac) conductance, we have 
measured dc conductance, because the cycles applied in 
the experiment [2^, 5.0 x 10~^ s, corresponding to 10^ 
steps of our MD are too long. We also carry out the MD 
simulations with A/i^ = 0.1, 0.4, 1.0, and 2.0 meV, and 
then obtain much the same conductance for A/i^ = 0.1 
and 0.4 meV as A^^ — 0.2 meV, and confirm no step- 
like structure in conductance for A/x^ = 1.0 and 2.0 meV 
(see Supplemental Material) For the /ii^ we set, the 
averaged density n over the channel is 2.32(±0.02) x 10® 
cm~^, where we enclose the standard deviation for vari- 
ation with Vg in the parentheses, and use its notation 
throughout this paper. Thereby, the 2D electrons be- 
long to a liquid state for the bulk system [5§| because 
the plasma parameter T = ^/irn / k^TK = 37.6 is much 
smaller than the critical plasma parameter (60j Fc = 137. 



E. Conductance calculation 

Under the above setup, we calculate the electric cur- 
rent, density distribution, potential energy, fluctuation 
of potential energy, etc. in a steady state after 500,000 
steps from the initial state (see Appendix |B] for details) . 
In order to obtain the conductance 



G{Vg) = -e/(yG)/AM, 



(9) 



as a function of Vq, we compute the electric current in 
the y direction in the channel at each Vq, 



I(Vg) 



2D 



(10) 



under the assumption — Afi^ . In this paper, (• • ■) 
represents both the ensemble and the time averages, 
where the ensemble is generated with the different ran- 
dom seeds for the initial velocity distribution. Through 
the current calculation, it is confirmed that the replace- 
ment of D in Eq. by j/min or ymax lowers / at each 
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Vg only by 1% or 2%. The contribution of to G is 
negligible because the corresponding quantity 

AA° - ( / - / ]d'r /WfcBTK In [Xln{r)] (11) 

iJy<0 Jy>0) 

is of the order of 10^'^A/i^. We also calculate G in the 
system with the reservoirs widened from the 20 x 20 fim^ 
squares to the 20 x 25 /im^ rectangles. G in the widened 
system is in the range of the error bar of G in the original 
system all over Vq. Moreover, we calculate G of the 
widened system with H = 20 fj,m, in which the obtained 
G in the widened system is in agreement with G in the 
original system. 

III. RESULTS 
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FIG. 3. (Color online) (a) The dc conductance G obtained 
from the MD simulations versus the gate voltage Vg, where 
each error bar represents the standard deviation of the en- 
semble average for each G. Here, Vg(1) = 0.15 V is the 
gate voltage at which the potential barrier Eb disappears, 
and G(1) is G at Vg(1). The solid line is the approximate 
conductance G in Eq. (|15[) . The vertical dashed lines and 
the number between them indicate the number of the stream 
lines organized by electrons in the gate, observed directly in 
the density distribution n(r). (b) dG/dVc with respect to 
Vg (red triangles) and the derivative of five-point unweighted 
smoothed G (red solid line). Here, the range AVs = 0.202 V 
between the vertical dashed lines runs from the first to the 
fourth peak in dG/dVa- 

In this section, we present the obtained results from 
our simulations. The results clarify the origin of the 
steps [26], and give clues for understanding the electron 
dynamics behind the steplike conduction. 



(«) ib) 




FIG. 4. (Color online) The spatial distribution of the effective 
potential energy 'l?(r) near the gate at Vg = (a) 0.09 V, (b) 
Vg = 0.18 V, (c) Vg = 0.25 V, and (d) Vg = 0.31 V. Here, 
the blue spheres are placed at an electron configuration of a 
step in our simulations. 
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FIG. 5. (Color online) The potential barrier Eb in Eq. (fT2|) 
(red triangles), and the standard deviations of temporally 
variational potential energy at the positions of the center of 
the gate, ctg (black circles), and in front of the gate, crp (blue 
rectangles). 



In Fig. [31^ a), we plot G under the time average over 
4, 300, 000 time steps and the ensemble average over 170 
different initial conditions, where G exhibits weak step- 
like increases. Moreover, dG/dV^ in Fig. [3]Jb) shows a 
characteristic oscillation similar to that in Ref. [26| . 

The insulation for low Vg is due to the existence of a 
potential barrier in the gate; on the other hand, the in- 
crease in G for high Vq is for expansion of the width of a 
constriction in the gate. Figure |4] depicts the spatial dis- 
tribution of the effective potential energy $(r) near the 
gate [see Eq. (|B2p ]. We should keep in mind that ^'(r) is 
the sum of the confining potential and averaged electro- 
static potential from self-organizing distributed electrons. 
As seen in Fig. (H^a), the existence of the high-energy 
barrier can be verified for low Vq. Figures [Hb)-(d) also 
display the decrease of the barrier height and the increase 
in the constriction width as Vg increases. 

With the aid of $(r), we introduce the energy barrier 
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FIG. 7. (Color online) The mean nearest-neighbor distance 
r2, r^, f4, and fs in the x direction among the electrons or- 
ganizing the two-, the three-, the four-, and the five-electron 
current, respectively. 



FIG. 6. (Color online) The density distribution n(r) near the 
gate at (a) Vb = 0.05 V, (b) Vq = 0.f5 V, (c) Va = 0.24 V 
and (d) Vb = 0.30 V. The hues observed in (b), (c), and (d) 
are formed by sequential single-electron flow. 



defined by 

i?B(^G) = ^(rG) - / d2r/(r;5F)<I>(r), (12) 

JSf 

which represents the potential energy difference between 
the center position in the gate, Tq — (0,0), and the re- 
gion in front of the gate, S'f = {{x, y)\yF < 2/ < 0}, with 
Uf = —1.7 /im. Here, ?/f is selected so as to include an 
average position of the electron just in front of the gate, 
and fir;SF) = n{r)/N{SF) with ^(^p) = Jg^d'^rn{r) 
is the normalized single-particle distribution function in 
S'f. Figure [5] shows that Eb becomes zero at Vq = 0.15 
V (= Vg(1)) and the barrier exists for Vq < Vg(1). 
However, a tiny current exists at Vq = 0.05 V in Fig. 
(JJa), although the electrons cannot get over the barrier 
{Eb = 2.02 meV) because of the small kinetic energy 
kgTK = 0.103 meV at Tk = 1 .2 K. The origin of the 
current will be discussed in Sec. IIV Al 

The steplike behavior in G is roughly understood by 
examining the density distribution function n(r) [see Eq. 
(|B1I) ]. FigurelHillustrates typical density patterns of elec- 
trons near the gate, corresponding to the s— th step of G. 
For low Vg, as seen in Fig. [DJa), the electrons seem to 
stay in front of the gate but flow slightly through the 
barrier. We can directly observe that the electrons flow 
through the gate in one line at the first step of G [Fig. 
ini^b)]), two fines at the second step [Fig. [H^c)] and three 
lines at the third step [Fig. IDJd)]. Therefore, it is con- 
firmed that the steplike increases are not originated from 
the conductance quantization [1], but the effects can be 
attributed to the increment in the number of stream lines 
of electron flow in the constriction, as suggested by Rees 
et al. [1^. However, the smooth steps in G cannot be 
understood only with the discrete increments in the num- 
ber Nq of electrons to pass simultaneously through the 



gate. 

As a possible mechanism for the smoothing, we can 
indicate the observed temporal variation of Nq in our 
simulations. To clarify this mechanism, we calculate the 
mean nearest-neighbor distance 



.(^g) - 



— r 
2 
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(13) 



among the electrons forming the m-electron current, 
which is the current when m electrons pass side by side 
through the constriction [see Eq. (|B5[) ]. Here, the region 
ScL = {(2^,?/)||y| < 0.05 /J,m } is set in order to mea- 
sure the separations in the conflned direction. It is also 
to be noted that the average of in Eq. ([TH)) is lim- 
ited to the middle distance of the nearest-neighbor and 
the second-neighbor distance. From fm for m = 2, 3, 
4, and 5 plotted in Fig. [7l we find that the two-, three-, 
four-, and five-electron currents appear only in more than 
Vg = 0.13, 0.25, 0.32, and 0.36 V, respectively. These Vq 
are lower than the Vg at which the m lines are observed 
in n(r) [see Fig[3ja)]. In addition, we find intermediate 
states in which the m-electron current with long separa- 
tions and the (m -I- l)-electron current with short sepa- 
rations coexist. Since G in the intermediate states takes 
between to- and (to -|- l)-electron current conductance, 
the effect can lead to the gradual change in G. These 
results also reveal that the separations between electrons 
to pass side by side through the gate have a fluctuation 
margin about 0.1 /im from the average interelectron sep- 
aration f = 21 \pKfL = 0.741 /im. A major factor in the 
temporally fluctuational separations will be specified in 
Sec. HV^^ 



IV. DISCUSSION 

In this section, we first discuss the origin of the current 
under the existence of the potential barrier and the mech- 
anism smoothing the steps in the conductance. Next, 
we compare the conductance in our system with the ob- 
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served conductance in the device [26|. Finally, the ob- 
tained conductance from the LD simulations is briefly 
compared with that from the MD simulations. 



A. Density-fluctuation-affected transport 
properties 

In this section, we propose two mechanisms of smooth- 
ing the steps in G. We attribute the current for Vg < 
Vg(1) to the intermittent disappearances of the barrier 
due to temporal fluctuation of the potential at both 
Vq and Sy- This mainly results from the following 
two reasons. First, the standard deviation of tempo- 
rally variational potential energy in the channel, a = 
<Prf{v] \y\ < D)a{r) = 0.74(±0.03) meV with Eq. 
(|B3p . is much larger than ksTK- Second, the velocity dis- 
tribution in the y direction near the gate for Vq < Vg(1) 
has the nearly zero mean value (~ 10~* — 10~^ meV 
in terms of the kinetic energy), and thus does not de- 
viate from the equilibrium distribution. Therefore, the 
conductor-insulator transition in our system is not caused 
bypinning and depinning of electrons at the constriction 

From the conservation of energy, it is clear that elec- 
trons can pass through the gate only when the potential 
energy exceeds the barrier. Therefore, G up to Vg (1) may 
be approximately represented as G sa G = G(1)P(1), 
where G(l) is the conductance at Vg(1) and P(l) is the 
probability that the potential energy of the front elec- 
tron exceeds the potential energy at re . From the direct 
calculations in terms of the MD (see Supplemental Ma- 
terial) |5^, we verify that the temporal change of (?I)(r,t) 
almost satisfies the normal distributions as follows: 



{S[e~<l>{r,t)]\{Sr),] 



(14) 



as in the case of the fluctuating electric field [57[. Fur- 
thermore, it is probable that the temporal changes of 
(/)(r, t) are uncorrelated between the two points of space. 
This is because the potential of the electron just in front 
of the barrier is fiuctuated mostly by two electrons in the 
rear [see Figs. lUJa) and[6lja)]; on the other hand, the po- 
tential at tq is fiuctuated mostly by the electron of the 
other side across the barrier. Hence, G can be estimated 



as 



G 



G(l) 
G 

G(l) erfc 



de 



de'- 









r (^--eb)"i 


exp 


2^ 


exp 


L J 



Eb 



where 



G — erfc 



4) 



Eb 



(15) 
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FIG. 8. (Color online) The spatial distribution of potential 
fluctuation a{r) (upper figure), and the density distribution 
n{r) (lower figure) at Vg ~ 0.05 V in the region of —10 fim 
< y < fj,m. 



is the probability for Eb ~ 0; crG(VG) — "'(re) and 
o'f(Vg) — Jg_^ d^r/(r; S'F)cr(r) are, respectively, the stan- 
dard deviation of the potential energy at tq and in 5f 
(see Fig. [S]); and erfc [a;] is the complementary error func- 
tion [solid line of Fig. Ela)]. Although G is estimated 
on the basis of only the one necessary condition, never- 
theless G is in good agreement with G calculated in our 
simulations. 

It should also be noted that the rising of G takes place 
at Vg at which Eb, ctf, and ctg satisfy the relation 



Eb ~2crF + 2crG, 



(16) 



-Eb=0 



where tails of the potential distribution at the two differ- 
ent points touch just each other. Therefore, this may give 
a method of an approximate estimation of the potential 
fluctuation from the potential energy. 

Concerning the choice of y-p in Sp [see Eq. p^ ]. the 
discussion on G holds well for the range of —2.3 /im 
^ Uf < —1.3 /zm in which the deviation |G — G| at each 
Vg is less than 0.04 MQ''^. Outside the range, G be- 
comes discrepant from G because of large spatial vari- 
ation in o'(r). Figure [8] illustrates that cr(r) increases 
in incommensurate regions in which two structures with 
the different number of the lines in n{r) are frustrated 

[il[ill6il,|6l. 

As seen in Sec. Illlj the temporal variation of A'^g arises 
from the fluctuational interelectron separations in the 
confined direction. The variational separations can be 
understood from an effect of the density fluctuation. The 
root-mean-square displacement Sx in the x direction es- 
timated from the harmonic approximation by equalizing 
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SlV^Mr) to kBTK/2 [5ll63|, i.e 



fee 7k 



2V20i(f)' 



(17) 



is given by Sx — 0.119 /xm, which is comparable with the 
fluctuation margin. The observed crystal-hke ordering in 
the confined direction in the density patterns in Fig. [5] 
also supports of our estimation. Therefore, we attribute 
grounds for the temporal change of Nq to the vibrational 
behavior. 

In the above discussion, we suppose that the smooth 
rising from the insulating state to the first step in G is 
due to the potential fluctuation, and the smooth steps 
in G are caused by the vibration in the confined direc- 
tion. These dynamics can be commonly attributed by 
the density fluctuations. Therefore, the magnitude of 
the density fluctuation seems to determine whether the 
steps in G can be observed. 



B. Comparison of our result with the observed 
conductance 




1.50 



0.25 



0.50 0.75 1.00 
Vg(V) 



1.50 



FIG. 9. (Color online) (a) The measured conductance Ge 
in the experiment (blue solid line) [2^ and G{Va) in Fig. 
Ha) transformed into G't(Vg) = G{[Vg + Vg(1)] AVs/AT/e - 
Vq (l))/(? (red triangles) vs the gate voltage Vq. Here, Vg(1) 
and AVs are given in Fig. [S] and AVe and Vq (1) are shown 
below, (b) dGE/dVo (blue solid line) ^ and dGx/dVo (red 
solid triangles), with respect to Vg. Here, the range AVe = 
0.750 V between the vertical dashed lines runs from the first 
to the fourth peak in dGE/dVa- Vq{1) = 0.63 V indicates 
Vg, corresponding to the experimental counterpart of Vg(1) 
in the simulations. 



In this section, we illustrate the qualitative consis- 
tency between our and experiment conductance in order 
to stress the suitability of our approach for the electron 
system in the device [26"]. 

The scaled and shifted G are similar to the measured 
conductance Ge in Ref. [26| as seen in Fig. [HI First, we 
select a scale factor g = 1.414 so that the magnitude of G 
at the first step fits into that of Ge. We confirm that the 
magnitude of G is dependent onjii or electron density 
as with the experimental results [20]. 

Second, we select a constant AVe/AVs = 3.3 as a scale 
factor for G in the Vg direction, where AVs and AT^e are 
the amount of increase in Vq with an increment of Nq 
from 1 to 4 in the simulation and the experiment (26| . 
respectively [see Figs, ^h) and IHb)]. The choice of 
the scale factor stands on the consensus that a steplike 
increase in conductance is determined by Nq. Because 
AVe decreases with decreasing electron density, as ob- 
served in Ref. [l^, and the density n = 2.32 x 10^ cm^^ 
in the simulations is lower than the experimental density 
fiE = 1.5 X 10^ cm-2, the result of AVs < AVe is 
reasonable. 

Third, we shift G so that Vg(1) agrees with the 
first-minimum-gate voltage Vq{1) which corresponds to 
the minimum between the first and the second peak in 
dGE/dVo [see Fig. [D^b)]. Because the guard voltage 
0.62 V in the device corresponds to V in our system, 
Vg ( 1 ) is practically larger than V^f ( 1 ) by { V^f ( 1 ) - . 6 2 } - 
{Vg(1) - 0} = 0.135 V. The shift is also valid. Indeed, 
the shift of conductance into the lower Vq direction is 
observed as electron density increases [26| . 

As a result of the above discussion, G transformed 
into Gt in Fig. [UJa) is almost in agreement with Ge 
for Vq > V(f (1), but is not good for Vq < V^{1). The 
disagreement for Vg < Vq{1) is also reasonable, for the 
amount of increase in Vq with the growth in conductance 
from the threshold of current flow to the first s tep is al- 
most invariant with respect to electron density [26|. For 
Vg < Vg{1), the shifted G without scaling of AVe/AVs 
is actually coincident with Ge for Vg < Vq'(I). 



C. Langevin dynamics simulation 

We also calculate the conductance in the molecular dy- 
namics under the Langevin thermostat, which is known 
as the Langevin dynamics, to check how the results de- 
pend on our choice of the Nose-Hoover thermostat. The 
detailed description on the method is given in Appendix 

o 

Figure [10] is the corresponding plot to Fig. ^a,) for 
the LD simulations. The obtained Gl in Fig. [TO] un- 
der the averages of 4, 300, 000 time steps and 16 different 
initial conditions reproduces the steplike conductance as 
observed in Ref. [26!|, and shows qualitatively similar 
behavior to G in Fig. [3]Ja). However, we can find the 
sharper rising to the first step of Gl than that of G. 
This deference reflects on the observed smaller potential 
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fluctuation in the LD simulations than that in the MD 
simulations (see Sec. IIV A|) . Because the characteristic 
time scale of the relaxation of electron motion is 1/^, 
which is introduced in Eq. (CI) 't^, the reduction in 
potential fluctuation may result from the damping. A 
more detailed comparison between MD and LD will be 
discussed in future work. 



range from the paper, which makes a study of the trans- 
port of the Wigner crystal [gl], the electron dynamics 
clarified by our investigations also provides information 
to help them understand the mechanism. 
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FIG. 10. (Color online) (a) The dc conductance Gl obtained 
from the LD simulations versus the gate voltage Vg, where 
each error bar represents the standard deviation of the en- 
semble average for each Gl. 



The authors thank D. Rees and K. Kono for fruitful 
discussions. This work was partially supported by the 
Ministry of Education, Culture, Science and Technol- 
ogy (MEXT), Japan (Grant No. 21540384), and by the 
Global COE program. The Next Generation of Physics, 
Spun from Universality and Emergence, from MEXT 
Japan. The author also thanks the Yukawa International 
Program for Quark-Hadron Sciences at Yukawa Institute 
for Theoretical Physics, Kyoto University. The numer- 
ical calculations were carried out on the SR16000 and 
Altix3700 BX2 at YITP in Kyoto University. 

Appendix A: Estimation of resistance due to the 
electron-helium gas atom and the electron-ripplon 
scatterings 



V. CONCLUSION 

In conclusion, we reproduce the steplike conductance 
in a classical 2D electron system confined in the point- 
contact device on liquid ^He by the molecular dynamics 
simulation. Conductance in two-dimensional electrons 
confined in the shape of a microchannel with a point 
contact by electric fields is calculated as a function of 
the confinement strength of the contact. It is confirmed 
that the number of the steps corresponds to the number 
of stream lines formed by the flow of self-organizing dis- 
tributed electrons at the contact. This result supports 
the expectation by Rees et al. [26] in which the conduc- 
tance does not originate from quantum effects but from 
classical many-particle effects at the contact. We ver- 
ify that a potential barrier exists in the contact for the 
stronger conflnement than that at the first step, and the 
barrier disappears for the weaker confinement in which 
the constriction width of the contact increases with weak- 
ening confinement. In the strong confinement, the rising 
of the conductance can be attributed to intermittent dis- 
appearances of the barrier due to temporal fiuctuation of 
the electrostatic potential. In the weak confinement, the 
number of electrons to pass simultaneously through the 
contact is incremented through the intermediate state in 
which the number varies temporally. Based on these re- 
sults, we suppose that the density fiuctuation smooths 
the steps in the conductance. 

Note added. Recently, a similar paper on the same 
subject based on both the experiment and MD has been 
published [64] , Although our work differs in temperature 



In this appendix, we investigate the transport process 
in the channel, as in Ref. [2^. It is demonstrated that 
the electron-helium gas atom and the electron-ripplon 
scatterings are ineffective in the experimental resistance 
i?E > 5 X 10~^ MSI for the system. 

For the bulk 2D systems, the effective collision fre- 
quencies, which are defined by the relaxation time in the 
Drude formula, are given in Ref. ^o[. For the setting 
in Ref. [11]: T = 1.2 K, zho = 1-5 ^m, Fj. = 6.5 x 10^ 
V/cm, the values of the frequencies are = 5.6 x 10^ 
s~^ for the the electron-helium gas atom scattering, and 

= 1.4 X 10^ s~^ for the electron-ripplon scattering un- 
der the single-electron approximation. We also confirm 
that the predominance of fg upon the resistance is un- 
changed even if the electron-electron scattering effect on 
v-c is considered under the complete control approxima- 
tion. Since these frequencies make the mean free path to 
be of the order of 0.1 /ini, the system is diffusive |ll| - [l7| . 

For the transport in diffusive regions, the macroscopic 
geometrical effects as in the Maxwell resistance are dom- 
inant [13]. Based on the Poisson equation and the Ohm 
law, the geometrical effect of the bottleneck [10 ^m 
< \y\ < 15 /xm; see Fig. ^a) and Ref. [2^] upon the 
resistance is calculated in Ref. [46| as 



TT 



* ^ 11 



tan 



tanh ^ 



log 



wt - wt 



(Al) 

where (P = m{i'g + v^) / n^e'^ is the bulk resistivity with 
the reported density ue = 1-5 x 10^ cm~^ in Ref. f26j . 
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and Wn = 10 and Ww = 2Wn are the widths at 
the narrow and the wide section of the bottleneck, re- 
spectively. Therefore, the resistance in the channel is 
given by i? = p^L/W^ + 2i?B = 4.7 x IQ-^ MTt, where 
L = 20 /im is the channel length without the bottleneck 
portion. Note that R is 10% of i?E. 

In addition, we take into account the boundary scatter- 
ing effects to the resistivity pj] . In the boundary scatter- 
ing theory [65l - l70j . the steady-state Boltzmann equation 
on the coordinate system in Fig. [H^a) under the free- 
electron approximation is given by 



dr 



(A2) 



which can be solved under an ideal hard-wall boundary 
condition, where ^'o is the equilibrium electron distribu- 
tion function to be the product of the constant density 
and the Maxwell-Boltzmann distribution function, VPi is 
a deviation of an electron distribution function from ^E'o, 
Fuis an in-plane driving field, and r is a relaxation time 
[6l,[7l|,[2l- In Eq. (jM]), we simply adopt r = l/(i/g + i/r) 
because the elastic scattering is dominant. We also set 
the boundary at a; = Wn/2 and x = —Wn/2 to sim- 
plify the argument. For the diffusive boundary scattering 
[65l [6^ , the resistivity p is given by 



2kBT 



T 



2fcBT n 

2 / aVxVxexp 

rmr Jq 



T 



1 



2kBT Vx 



(A3) 



This result leads to p = l.OSp", and thus, the effects 
little affect the channel resistance. Therefore, the corre- 
lation effects between electrons in the external potential 
are dominant for reproducing the value of -Re- 

I 



Appendix B: Calculation method of physical 
quantities 

The one-particle quantities we calculate are the density 
distribution function averaged over the 0.1 pia square 6r 
centered around the position r: 

n(r) = ^^5(r-r,(t))^, (Bl) 

the average potential energy of electrons in dr, 

$(r) = (^(r„t)|(<5r)i), (B2) 
and the standard deviation of temporal fluctuation of 



{6r),)}\ (B3) 



a(r) = {(|</.(r„i)-ci>[r,(t)] 

where (• • • |((5r)i) represents the ensemble average under 
a conditional time average for electrons in Sr, when elec- 
trons are present in Sr. Here, Sr is created by dividing the 
xy plane into 0.1 pm square mesh, and n in ((5r)„ denotes 
that the averaged quantity is the n-particle quantity. The 
conditional time average is defined by 



(B4) 



-4— dt f <fr' y (5 [r' - r,(i)] • • 

where is the calculation time, and 

<p(Jr) = /*" dt f dh" y S [r" - rk{t)] 

Jo JSr ^ 

is the total presence time of electrons in Sr. We also 
confirm that the quantities are almost in agreement with 
those calculated with 6r to be 0.05 pm. squares. In ad- 
dition, we use another conditional average for a two- 
particle quantity as follows: 



(•••|(^)™) = (^^;(|^^*"dt'5(m-iVsW) l^d'r' J^d'r"J26[r' -r,{t)]S[r" ^r,{t)]---\ , (B5) 

where {• • •)e is the ensemble average, Ns{t) is the number of electrons in an area 5* at time t: 

Ns{t) = / d2r'V^[r'-rfc(t)], 

and t2p{S; m) is the product of the number of pairs among m electrons and the presence time when just m electrons 
are present in S, 

t2p(5;m)= r dt6[m-Ns{t)] f dV / dV V ,5 [r' - rfe(t)] <5 [r" - (i)] . 

Jo JS JS 



I 

Thus, (• • • 1(5*)™) describes the ensemble and the condi- m electrons are present in S. Incidentally, all of the in- 
tional time averages, among m electrons in 5, when just 
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troduced quantities in this section are time and ensemble 
averaged over 1,200,000 time steps and 32 different initial 
conditions, respectively. 



Appendix C: Langevin dynamics simulation method 

In this appendix, we briefly explain the method of LD 
simulation to calculate the conductance of electrons [z^j ■ 
The equation of motion for the electrons in the LD sim- 
ulations is given by 



(CI) 



associated with Eqs. ([7]) and ([8]), where consists of the 
first, third, and fourth terms of the right-hand side in Eq. 
([5]), ^ is the friction constant, and F^ is a random force, 
reproducing the thermal noise with zero mean value and 

2mkBT^dijSapS{t — t'). If we 



variance {F^{t)F!l {t')) 



regard ^ as the collision frequency through the Einstein' 
relation, we may estimate £^ as £^ = + = 7.0 x 10^ 
s^^ (see Appendix VK\ for fg and ^-r). For the value to be 
^ <C 1/ts, we have no reason to neglect the inertial term, 
where 1 /ts is the characteristic frequency for the correla- 
tion between electrons. To integrate Eq. (jCll) . we adopt 
the Ermak's approach tH, 3, which performs properly 
the stochastic integration [7a]. This LD simulation is re- 
duced to an overdamped LD simulation for t^S^ 3> 1 , and 
our MD simulation without the Nose-Hoover thermostat 
in the limit ^ — >■ 0. The other numerical methods and 
the simulation setup are unchanged from those of the 
MD simulations (see Sec. HI]) . 
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